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SECTION 1 « INTRODUCTION 

The interactions of dynamics, photochemistry and 
radiation in the stratosphere can he described^ by the continuity, 
momentum and energy equations. A more comprehensive model for 
stratospheric transport theory has been developed for the purpose 
of aiding predictions of changes in the stratospheric ozone content 
as a consequence of natural and anthropogenic processes. This 
model is time dependent and iv dependent variables are zonal 
means of the relevant meteorological quantities which are 
functions of latitude and height. The report sets out the detailed 
formulation of a numerical model both in physics and mathematics. 

A set of fundamental dynamical equations will be described in 
Section 2. Section 3 describes the numerical method used in the 
integration. The results and discussions are presented in 
Section 4. 
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SECTION 2 - BASIC EQUATIONS 


The model is an integration of an extensive set of 
fundamental equations which can be used to predict the future 
state of the atmospheric circulation from knowledge of its 
present state. In the log-pressure coordinate system (x, y#^) 
in which x » distance eastwards, y = distance northwards, and 
^ = log (Pq/p) ; with p^ = 1000 mb and p = pressure, u, v, and w 
are the three corresponding velocity components respectively. 

In this log-pressure system the eastward component of the momentum 

^ . .2 

equation is: 

u + u u + V u + w u = tan<^ + v sin^- ^ 

dy a dx 

where a = radius of the earth, = rate of rotation of the earth, 

(p = latitude, Z is geopotential. The thermodynamic equation is: 

where 0 = T (p^/p) , T = temperature, q is the rate of change of 
potential temperature, k is the ratio of gas constant for dry air 
R to specific heat of air at constant pressure Cp. The equation 
for conservation of matter is: 


du . dv 

“57 


V 

— ^ tan (p + 



w 


0 . 


(3) 
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If Xi is the mixing ration of a substance in the 


stratosphere, then the eq. for the time change of may be 
written in eq. (4) where Pj|^ is the production rate which includes 
the production and loss terms. 



9X1 

3x 


+ V 


3Xj 

3y 


+ w 


1. = P 


(4) 
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Our aim is to deal explicitly with zonal mean 
quantities only. Accordingly, let ai over bar denote the zonal mean 
and a dash denote the departure therefrom; then equations (1) to 
(4) yield, after taking an average around a latitude circle:^ 


d 



COS(p) 



+ 


^(Wf) 


H 


(5) 


-“ ( 6 e ^cosQ?) 


-J^(V e ) 


d 


(W ^ ) = q e”^ cos ^ + Q 


(6) 


3y 


V + 



0 


(/) 


-y^(Xi e"*^ cos<p) + -^(V Xi) + Xi) =Pi6 ^cos<p+ F (8) 

where T = (SI a cos (/> + u) a coscp 

„ = - a., 77) - -32-(~) 

Q = - 3^(771 - -3|-t77 ) 

F = - 37<V'Xi) - 3^(w'Xi') 

V - V e ^ cos (p 
VJ = w e ^ cos(p 

Note that the continuity equation (3) takes the simnle form (6) which 
is the advantage of using the pressure coordinates. 
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I 

Equations (5) to (8) may become a complete set by 
adding the y component momentum and hydrostatic equations. The 
hydrostatic approximation has been used in the ^ direction. For our 
own purpose of not introducing the gravity waves we introduce 
an approximation to the y component momentum equation that the 
zonal mean wind is in geostrophic equilibrium with the zonal 
mean geopotential field. Combining this with the hydrostatic 
equation leads to the thermal wind relationship: 

_Jl_ r - - R exp(-kp .5. 

3^ 2 ft tan<^ ay 

Equations (5) to (9) can be regarded as a closed set 
which completely determines the relationships among the dependent 
variables t, V, W, Xi. The other variables can be written 
in terms of the relevant quantities and are held as functions of 
time, latitude and altitude. 
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SECTION 3 - NUMERICAL METHOD 


It is convenient to make the equations dimensionless. 
The scaling quantities are given a subscript s. Equations (5) to 
(9) become: 


— (re cos^) + (TV) + (t W) 


3t 


H 


( 10 ) 


{$ e cQ'^'p) + 
3 w 


-i- f = - e 
H 




tan^ 


Tv <® ''’ * iV ‘® 

^(p 


- 

q e ^ cos g>+ Q 


(12) 


V =- 




w = 

3^ 


(13) 


^ (X e cose?) + -5— (V X) + (W X) = P e ^®^cos</?+ F 

d(p at ^ 


3t 


where H and Q are the divergence of edciy momentum and eddy heat 
flux, respectively. 'I' is a stream function which is introduced 
to rewrite equation (7) in the form of equation (13) . The method 
of solution is to eliminate the time derivatives from (10) and 
(11) using (12). Substituting for V an'd W from (13) gives an 
equation for 4'. 

0(pe 




6tc , _ 0(pe ^'^^s(l + sin^c?) , 

t ~^cp<p~ rep^ -2 ^ ^ 

V' L sin <p 







+ 
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je ^^scos^ + Q tan^Jj /tan^ (15) 

The 4 and 99 used as subscripts denote differentiation with respect 
to that quantity. It nay be solved with suitable boundary 
conditions for ^ which may be used in ( 10 ) and ( 11 ) to get new 
values of t ,$ at the advanced time step. Equation (15) is then 
solved at the new time value and the process repeated. 

4 

A Crank-Nicholson implicit finite difference scheme 

is used in the numerical calculations. The Alternate Direction 
. . 4 

Implicit method is adopted so that the finite difference equation 
can be solved by simply inverting a tridiagonal matrix in each 
of the two dimensions sequentially. The nonlinear coefficients 
in the equations can be taken into account by using a predictor - 
corrector procedure. A more detailed description follows 

Eq. (15) is an elliptic type, second order partial 
differential equation which can be solved by running a marching 
method for the time^dependent equation to a sufficiently large 
time. Equations (10), (11) and (14) are time-dependent first 
order partial differential equations. A general form for the set 
of equations can be written as: 

= A + B + C + D + E + F (16) 

The ADI scheme for equation (16) is the following: 

+ At [^A6| (S^'^'^) + B C50(s") 

+ D5^ (S”) + f1 (17) 


3-2 



gn+1 ,gn ^ g^2^gn+>jj g gn^ + C (S^'*’^) 

+ D6^ (S""^**) + E5<^ (S""*"^) + pj (18) 

where ; 

■®| ^jk ” ‘®j+l,k" * ®J,k * ®j-l,k*''^5^ 

®jk “ '®j,k+l" ^ ®j,k ■'^ ®j,k-l’^‘^'^ 

«{ ®3k - - S^.^,)/2A{ 

«<. ®jk “ <®j,k*l - ®j,k-l’^"*''’ 

A simple tridiagonal matrix inversion will roIvo equations (17) 
and (18) to give the solution of equation (16). 


j-3 


Eq. (15) is an elliptic type, boundary value problem 
so that the crossing termi//^^ is not a dominant term which can 
be approximated at the old time step as shown in eqs. (17) and 
(1C). This enables us to use the Crank-Nicholson scheme and 
the ADI method to solve eq. (15) without rotating the coordinates 
to get rid of the crossing term. We have ta)cen the average for 
each two consecutive values of in the process of getting fast 
convergence. In practice Eq. (15) has been solved with less 
than 50 total time steps to reach the steady state solution. 

The values of ij/ are then used to determine V and W 
which are then substituted in Eqs. (10) , (11) and (13) for solving 
T, e and X at the advanced time step. This process has been repeated 
fcr the time dependent solution. 
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SECTION 4 RESULTS AND DISCUSSIONS 


The terms in the riaht hand side of Eq. (15) are source 
functions which drive the wind fields. Q is the sxunination of 
the rate of heating due to absorption of solar radiation and the 
rate of cooling by thermal radiation. Eddy fluxes of heat and 
momentum are treated by large-scale diffusion coefficients. For 
these fixed driving terms the time-dependent solutions of the 
set of equations depend strongly on the initial data available. 
Therefore we have to find the steady state solution first before 
seeking any time -dependent solutions. 

For the winter northern hemisphere case the generated 
zonal (u) , meridional (v) and vertical (w) winds are shown in 
Figure (1) , (2) and (3) respectively. The agreement between the 
observed and calculated data may be considered fully satisfactory. 
Figure (4) shows the temperature distribution from south pole 
to north pole (x axis) as well as from ground to about 55 km 
(Y axis) approximately. As we mentioned early H, Q and q are 
the three source functions which are shown in Figure (5) to 
Figure (7) respectively. From Eq. (15) it is possible to regard 
^ as the sum of separate contributions by a different term on 
the right hand side. Thus one may set to zero all the terms on 
the right-hand side of (15) except those in q and the solution 
then gives that part of the mean circulation induced by the net 
heating. The meridional and vertical winds arising from heating 
only are shown in Figure (8) and (9). It shows the heating plays 
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an important role in the stratosphere but not in the troposphere. 
Figure (10) and (11) show the results arising from eddy heat 
fluxes. Finally Figure (12) and (13) are due to eddy momentum 
fluxes alone. The structure of the wind fields in the troposphere 
are heavily determined by eddy fluxes of heat and momentum. 
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Figure Captions 


For the figures shown below, a northern hemisphere 
winter case (December to February) has been chosen. The x 
coordinate shows the latitude from south pole to north pole. 

The y coordinate gives the pressure and altitude in the follow- 
ing way. 
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Figure 1. Zonal Wind (tn/sec) 

Figure 2. Meridional Wind (m/sec) 

Figure 3. Vertical Wind (m/sec) 

Figure 4. Temperature ( k) 

Figure 5. q (net heating) (dimensionless) 

Figure 6. Q (eddy heat flux) (dimensionless) 

Figure 7. H (eddy momentum flux) (dimensionless) 

Figure 8. Meridional Wind due to q only 

Figure 9. Vertical wind due to q only 

Figure 10. Meridional wind arising from Q only 

Figure 11. Vertical wind arising from Q only 

Figure 12. Meridional wind due to H only 

Figure 13. Vertical wind due to H only 
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Figure 6 
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SECTION 1 PHYSICS AND MATHEMATICS 


The interactions of dyneunics, photochemistry and 
radiation in the stratosphere can be described^by the continuity, 
momentum and energy equations. In the log-pressure coordinate 
system the non-dimensional forms of these equations become 


3 

3t 

(T 

e ^ COS0) + 

3 

30 

{tv) + 

3 

(TW) = H 

(1-1) 

3 

3t 

(0 

e ^ COS0) + 

_3_ 

30 

(0V) + 

3 

H 

(0W) = q e ^ COS0+ Q 

(1-2) 


H 


T = - 


tarup d<p 


+ _? = o 

34 


— (X COS0) + (VX) + ^ 


(1-3) 

(1-4) 


- 4 * 


(WX) = P e ^ cos<p + F (1-5) 
2 


30 ' ' 34 

where the notation of Harv;ood and Pyle'' has been used. This is 
a complete set of equations for the unknowns v, w and x. 

Equation (1-4) enables us to introduce a stream function such 
that 


V 




( 1 - 6 ) 


The method of solution is to eliminate the time 
derivatives from (1-1) and (1-2) using (1-3). Substituting V and 
V'J from (1-6) gives a equation for 0 


-T0 + 2 0^0 



( tan0) 


-1 
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+ e"^^* (1 + sin^<p) 


(sin<f>) ^ S(f> ^(f>^ 


+ k( g e“^«* (tan</»)"^ 




~ + e“^(* (tan</>) 


-1 


-e* 


e COS0 q^ + P<^> + Q tan^ 


(1-7) 


For suitable boundary conditions Eq. (1-7) can be solved for 

which gives us V and W. Then (l-l) , (1-2) and (1-5) are solved 

for T,.^andxat the advanced time step. (1-7) is again solved 
N 

at the new time value and the process repeated until the solution 
become stationary. 

A general form for the set of equations can be 

written as 

= A + B + C S00+ D + E S0 +F (1-8) 

A Crank-Nicholson implicit finite difference scheme^ is used in 
the numerical calculations. The Alternate Direction Implicit 
method is adopted for this two dimensional problem so that 
the finite difference equation can be solved by inverting a 
tridiagonal matrix in each of the two dimensions sequentially. 

The ADI scheme for (1-8) is the following 


+ At 


A 6^ (S"'^^) + B S^0 + C 60 (s"") 


+ D 6^ + E 50 (S*^) + F 


+ At 


A 6| + B S^0 + C ^ (S'^'^^) 


+ D 6^ (S'^'^^-') + E 60 (S“’") +• Fj 


n+1, 
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where the following notations are used 


4 ®j,k “ ’®j,k+l ■ ^ ®j,k ®j,k-l’/ 
®? ®j,k ' ‘®j+l,k ■ ®j-l,k'/2A| 

% ®j,k ” '®j,k+l ■ 
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SECTION 2 PROGRAM OUTLINES 


2.1 Main 

The developed GSFC 2D model (latitude, altitude) 
introduces a better understanding of the coupling between 
stratospheric chemistry and dynamics. The program is constructed 
in American National Standard FORTRAN on a large IBM S/360-91 
compute: with a region less than 400K eight-bit bytes. It consists 
of main program, seven (7) subroutines and data. In the main 
program, the physical constants are defined and the mesh sizes 
in each of the two dimensions are chosen. The eddy diffusion 
coefficients, the initial profiles for zonal wind, temperature 
and ozone distribution are read in from the punched - card data. 

The five major subroutines accomplish the following tasks: 

(1) Subroutine HEAT deals with the calculation of heating and 
cooling rates of the atmosphere (2) Subroutine EDDY calculates 
the eddy diffusion terms due to heat fluxes and momentum fluxes. 

(3) Subroutine \^'JS computes the solutions for meridional and 
vertical winds. (4) Subroutine TVS solves the time-dependent 
solutions for zonal wind and temperature and finally, (5) Subroutine 
ANS prints out the solutions in two-dimensional format. A 
schematic diagram including the JCL cards is shown in the following 

///ZMHCCJ28 JOB { L4004 2 19 5A , T , AAOO'0 1 , 00100 1 ) , 620 , MSGLEVEL=1 
//*B11 

// EXEC FORTRS.NH,PARiM= 'OPT=2 ,XREF' ,REGION=400K 
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Main Program 
Subroutine EDDY 
Subroutine VWS 
Subroutine HEAT 
Subroutine TUS 
BLOCK DATA 
Subroutine TRIDJ 
Subroutine ANS 

// EXEC LOADER, REGION=400K,PARM='SIZE=400K' 

//FT05F001 DD DSN=K3 , L6CJM. ZMHCC . SI , DISP=SHR 
// DD *,DCB=BLKSIZE=2000 

/* 

Every subroutine will be described in great detail in the later 
sections. 

The M.K.S. System has been adopted in the calculations. 
Some of the important physical constants and parameters are listed 
and defined below. 

PI : 7T 

RA : radius of earth (m) 

2 

GRAV ; acceleration of gravity (m/sec ) 

XM : mass of air (Kg) 

BOLTZ ; Boltzmann constant, (k) (*^^^^/deq (k) ) 

XC : speed of light (c) (m/sec) 

XH : Planck constant (h) (dyne-sec) 

XK : ratio of gas constant to specific heat at constant 

pressure. ();) 

JJ : number of grid levels in vertical direction 


1 
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R 

CS 

t 

v;s 

KK 

KT 

DZ 

DY 

Y(J) 

EX(J) 

EXK(J) 

P(J) 

Z(J) 

SI 

CO 

TA 

DT 


angular velocity of earth 

^ 8 
Vs 

«s 

total number of intervals in latitude 
total grid points in 12 hour average 
grid size in latitude 
grid size in altitude 

4 

exp ( - s) 
exp (- kHs) 
pressure (mb) 
latitude (degree) 

t 

sine ( 
cosine (6) 
tangent (0) 
time step (sec) 


The following parameters are read in from the machine 
readable data stored on a CRBE file on disk: 

YYK : horizontal eddy diffusion coefficients 

YZK : off-diagonal eddy diffusion coefficients 

ZZK : vertical eddy diffusion coefficients 

TEMP : temperature (°k) 

ZA : wavelength (m) 
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-2 -1 

ZB : solar flux (m sec ) 

2 

ZC : absorption cross section of 0^ (m ) 

2 

zy : absorption cross section of Oj (m ) 

-3 

P03 : ozone density (in ) 

UU : zonal wind velocity (m/sec) 

It is then ready to call subroutines Heat, Eddy, Vws, 
TUS sequentially to solve the problem. 


2.2 Subroutine HEAT; Heating and cooling rates are 

calculated mainly from the 15 fim CO 2 band and the 9.6 /im band. 
The parameters used are defined below: 

RD : unit from radians to degrees 

C2 ; CH * XK *86400 *1000/i.oi3 x 10^ 


AL 


AiM 

P02 

0 

QQ 

XJ2 

XJ3 

Ql 

Q2 

O"'- 


XB 

XX 

XF 


g m/ k t 

density of air (m 
density of ©2 (m 
mixing ratio of 
heating rate (°k/day) 

J coefficient for ©2 
J coefficient for 

heating rate at different hour of the day 

XJ2 at different hour of the day 

XJ3 at dif fc-rcn h hour of the day 

column density for O 2 

column density for 

distance from the center of the earth 
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KT 


total number of points in 12 hours 


S4 

TAU 


: Ch (X) 

• r 0.(X) y*® N. ( 2 ) CH(x) 

The rate of cooling is written 


= XK2C02 - XK1C02 * 

XK203 - XK103 * 3^ (T) 


where B{T) 

27Tk-'T-' 

x-^ 

ch^ 

e'^'-l 

and X = k T 


S3 : 

cooling rate 

due to CO 

S6 : 

cooling rate 

due to 0^ 

SO : 

net heating 

rate 
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2.3 


Subroutine EODY; H and Q are the divergence of eddy 
momentum and eddy heat fluxes respectively. They are defined 

as 


H « 

30 


_ 3 _ 

30 


34 


C J: iL + K - ^ 

yy aus 30 yz 03 34 


e ' COS0 


K, 


^ 2 -i. + K — - — (_ii_)_jLL_] 

- 30 ^2Z ' 3z ^ 34 J 


yz a Wg 32 


YYK + YZK , 
30 34 J 


Z YK + 

30 


ZZK 


H J 


~ t* 

e ’ COS0 


e~^ COS0 


e ^ COS 


0 - 


30 


YYK 4^- + YZK -1.^—1 e'^^cos 


30 




_L 


+ I ZYK -4^— + ZZK - 4—1 e cos 


30 34 

where T is reoresentcd in the program by TH which is the average, 
between old value T and new value TT. Similarly 6 is represented 
by UH. The eddy flux convergence of mixing ratio is represented 

by PF and PL is the production rate. XK2 and XK4 are rate constants 
for the basic reactions 


O + O 2 + ^5 — ► 0^ + M 

and 0 + 0^ — — ► 2 0^ respectively. EN is the enhancement for XK4 in 
order to take into accounL the ozone sink due to other species in 
the atmosphere. 
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Subroutine WS : This subroutine solves a second order 

partial differential equation for stream function i// in the form 
A + B + C + D \l/^+ E V^+P = o 

where A = ' T 0 tan0 
B = 2 tan0 

D = "{ |gT 0 + e“^^* (— ^ + 1 ) tan 0 

' sin 0 

^ ^ ^ ^ ® ^ ^ ^s ^*0 

F =-(H^ + 4^H) tan0 - e“^^* (e"^*cos0 q 0 

+ Q 0 + Q tan0) 


The alternate-direction-implicit method' has been used to solve 
tnis two-dimensional problem. The finite difference equation for 
the first step is 




AA. 

3 

= -A + D 



BB . 

= -4- + 2 

A 


3 

At 



cc^ 

II 

1 

> 

1 

o 



YY . 

" s. . , 

+ 

1 

< 
1 

s „ r j_ 

3 

3-1, k 

L J 



+ S . , . 

3 ,k+l 

[c + eJ + 

"j,k-i [c - 


where 


B 

Q — A A I 


4A? A 0 rj+l,k-l ■ ®j+l,k-l '^j-l.k+X * ®j-l,k-l) 


The second step is 


* ST^ ^ + CC,. * ST, = YY 


j 


j fk+1 
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where AAj^ = (~ C + E)/2 


= sr 


+ C 


CCfc - (- c 

- SS._ 

- E)/2 

l,k h- 

d| ^ ss j 

» A 

1-2*1 + 

+ S • 

i A - n + S . . T-i — 

2A - C + 

D-l.k 

1 * 

D»K lAt 

J 


■ C El 


C ^ E ■ 



®j,k+l 

’T' ^ ~T~ 


+ F 


+ g 


Subroutine TRIDJ is called to invert the tridiagonal 
matrix. II is the total number of time step required to satisfy 
the criterion (S (J,K) - ST(J,K) ) / S (J,K) s .005 for steady 
state solution.^ Integer I is designed to average the solutions 
for every two time steps which speeds up the rate of convergence. 
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2.5 Subroutine TVS: Lquations (1-1), (1-2) and (1-5) are 

solved at the advanced tiin'j ctep. They have the general form 


S. e“^*cosd) = V + W + RHS 


b4> 

The finite difference equation for the first step in ADI scheme is 

AAu * SS . , . + BB. * SS. , + CC. * SS . , = YY . 

where AA^ = - W/4A^ 

BB^ = e ^ COS0 / At 
CCj = W/4A4 

YYj = (W/4AO + (e“^*cos0/A t) 

"^j+l,k (W/4A4) - (®j,k+l ■ ^j,k-l) 

The second step is 

"‘k * ^ * ®^j,k " * s-^j.k.r ’'’'k 

where AAj^ = - V/4 A<p 

-t* 

= e ^ cos 6 / A t 
CC^ = V/4A0 


4- RHS 


j + l^k"^ ^j+l,k 


yv = _y_ / \ _ ^ f 

‘k 4A0 \ j,k+l j,k-l / 4A^ y 

- ®=3-l,k - ®j-l,k) ^ e-«‘cos^./At + RHS 
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SECTION 3 COMMON 


In the program a double precision has been chosen 
for all real variables. 

A . Common /TU/ 


FORTRAN Name 

T(28,19) 

TS(28,19) 

TT(28,19) 

US(28,19) 

UT(28, 19) 
ii(19) 

Zl(28) 

Y(28) 

DZ4 

DT 

DTO 


Description and Units 

(absolute angular momentum) 

T (intermediate values beteeen 
s 

"■t Vl) 

Tt+1 


6 (intermediate values between 
s 

®t ®t+i 
0t+l 


latitude (degree) 
yy for O^ 

( - log (p*/1000 mb) ) 
4 * DZ 


real time step (sec) 
1/DT 


B. Common/HEA/ 


03(28,19) 


P(28) 

AB(19) 

PP(19) 


0 (intermediate 
s 

°t °t+4t 

pressure in mb 

not used 

not used 


values between 


QQ(2S,19) 


heating rate (°k/day) 



ZA(162) 

wavelength (m) 

ZB(162) 

-2 -1 

flux of solar photons (m sec ) 

AC{162) 

absoprtion cross section of 
0^ (m^) 

ZY(162) 

eUasorption cross section of ©2 (m 
-6 

X03 

CH/k/9.6 X 10 0^ band) 

XD 

27T / (c * 100. * h^) 

XC02 

CH/k/1.5 X 10 ^ (15 fim. CO 2 band) 

CH 

C * h 

BOLTZ 

k (Boltzmann constant) (dyne/®k) 

RA 

radius of earth (m) 

SD 

sine of declination angle 

CD 

cosine of declination angle 

PI 

7T 

GMK 

g m / k 

KT 

total points for 12-hour average 

KTT 

2 (KT-1) 

KTMl 

(KT-1) 

C . Conmon/W'?/ 

q(28,19) 

defined in (2.4) 

ZZ (28) 

YY for r 

SQ(28, 19) 

q (net heating rate) 

EXK(28) 

exp (- k 4 ^g) 

SI (19) 

sine (latitude) 

TA(19) 

tangent (latitude) 
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AA(162) 

BB(162) 

CC(162) 

YY(162) 

W(28,19) 

XX(28) 

ST(28,19) 

SH(28,19) 

U(28,19) 

5(28,19) 

55(28,19) 


AA . or AA, 

BB . or BB. 

3 ^ 

CC . or CC, 

J ^ 

YYj or YYj^ 

vertical wind velocity (m/sec) 
solution for subroutine TRIDJ 
solution for ’A' at t + At 
average of 'A' between t and t + At 
potential temperature 9 
solution of i// at t 

4's (intermediate values between *A'+ 
and P t+At) 


111 

HIT 

1111 

llllT 

DY4 

DZDZ 

DYDZ4 

DYDY2 

XK 

CS 

D. Common/EDD/ 

A(28, 19) 

B(28 , 19) 

C(2B, 19) 


loop for integration 
end point for loop 111 
loop for integration 
end point for loop 1111 
4 * DY 
DZ * DZ 
4 * DY * DZ 
2 * DY * DY 
1/Cp 


coefficient 

in 

eq. 

(1-8) 

coefficient 

in 

eq. 

(1-8) 

coefficient 

in 

eq. 

(1-8) 
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r. 

] 

D(28.19) 

E(28,19) 

YYK(28,19) 

YZK(28,19) 

•ZZK(28,19) 

PL(28,19) 

ZYK(28,19) 

W4(23,19) 

P02(28,19) 

P03(28,19) 

TEMP (28,19) 

AL(28,19) 

XJ2(28,19) 

XJ3(28,19) 

UU(28,19) 

F(28,19) 

Q(23,19) 

H(28, 19) 

UH(28,19) 

TH(23,19) 

0(28,19) 

PF(28,19) 

CO (19) 

EX (28) 

KKPl 

JJPl 


coefficient in eq. (1-8) 
coefficient in eq. (1-8) 
eddy diffusion coefficient 
eddy diffusion coefficient 
eddy diffusion coefficient 
production rate 
eddy diffusion coefficient 
density of air (m ^) 
density of O 2 (m ^) 
density of 0^ (m ^) 
temperature ( ° K) 

1/H (scale height) 

J coefficient for O 2 
J coefficient for O^ 
zonal wind velocity (m/sec) 
coefficient in eq. (1-8) 

Q in (H ft P) (dimensionless) 

H in (H ft P) (dimensionless) 
average of B between t and t + At 
average of r between t and t + At 
mixing ratio of 0^ 
eddy diffusion term for 
cosine (latitude) 
exp (-^*) 

KK + 1 
JJ + 1 


K 


yy 


yz 


zz 


zy 


I 
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DZ2 

2 * DZ 

DY2 

2 * DZ 

DY2 

2 * DY 

DZ 

mesh size for latitude 

DY 

mesh size for altitude 

KK 

total intervals in latitude 

JJ 

total intervals in altitude 

AR 

Q * radius of earth 
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List of Symbols 


T 


9 

H 

Q 

J 

V 


w 


a 

o 

4 . 

q 

t 

R 

u 

V 

w 

Q 

'<> 


absolute angular momentum 
= (Ua co3<t>+ u) a cos<^ 
potential temperature 

= - - iV ‘ZE’ 

= - - if 

Photodissociation rate 
VCOS0 exp (- ^ ) 
wcos(^ exp (- ^ ) 
radius of the earth 
pressure 

rate of change of potential temperature 
time 

gas contstant for dry air 
eastward veolcity component 
north’./ard velocity component 
vertical velocity 
rate of rotation of earth 
latitude 

stroain function for mean meridic^nal circulation 
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